/*
 * Copyright (c) <2008 - 2020>, University of Washington, Simon Fraser University, Bilkent University
 * All rights reserved.
 *
 * Redistribution and use in source and binary forms, with or without modification,
 * are permitted provided that the following conditions are met:
 *
 * Redistributions of source code must retain the above copyright notice, this list
 * of conditions and the following disclaimer.
 * - Redistributions in binary form must reproduce the above copyright notice, this
 *   list of conditions and the following disclaimer in the documentation and/or other
 *   materials provided with the distribution.
 * - Neither the names of the University of Washington, Simon Fraser University, 
 *   nor the names of its contributors may be
 *   used to endorse or promote products derived from this software without specific
 *   prior written permission.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
 * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
 * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
 * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
 * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
 * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
 * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 */

/*
 Authors:
 Farhad Hormozdiari
 Faraz Hach
 Can Alkan
 Emails:
 farhadh AT uw DOT edu
 fhach AT cs DOT sfu DOT ca
 calkan AT cs DOT bilkent DOT edu DOT tr
 */

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include "Common.h"
#include "CommandLineParser.h"
#include "Reads.h"
#include "Output.h"
#include "HashTable.h"
#include "MrFAST.h"

char *versionNumber = "2.5"; // Current Version
unsigned char seqFastq;

int main(int argc, char *argv[]) {
	if (!parseCommandLine(argc, argv))
		return 1;

	configHashTable();
	/****************************************************
	 * INDEXING
	 ***************************************************/
	if (indexingMode) {
		/********************************
		 * Regular Mode
		 ********************************/
		configHashTable();
		generateHashTable(fileName[0], fileName[1]);
	}
	/****************************************************
	 * SEARCHING
	 ***************************************************/
	else {
		Read *seqList;
		unsigned int seqListSize;
		int fc;
		int samplingLocsSize;
		int *samplingLocs;
		double totalLoadingTime = 0;
		double totalMappingTime = 0;
		double startTime;
		double loadingTime;
		double mappingTime;
		double lstartTime;
		double ppTime = 0.0;
		double tmpTime;
		;
		char *prevGen = getMem(CONTIG_NAME_SIZE);
		prevGen[0] = '\0';
		char *curGen;
		int flag;
		double maxMem = 0;
		char fname1[FILE_NAME_LENGTH];
		char fname2[FILE_NAME_LENGTH];
		char fname3[FILE_NAME_LENGTH];
		char fname4[FILE_NAME_LENGTH];
		char fname5[FILE_NAME_LENGTH];

		char outputFileName[FILE_NAME_LENGTH];
		// Loading Sequences & Sampling Locations
		startTime = getTime();
		if (!readAllReads(seqFile1, seqFile2, seqCompressed, &seqFastq,
				pairedEndMode, &seqList, &seqListSize)) {
			return 1;
		}

		loadSamplingLocations(&samplingLocs, &samplingLocsSize);
		totalLoadingTime += getTime() - startTime;

		if (pairedEndMode) {
			//Switching to Inferred Size
			minPairEndedDistance = minPairEndedDistance - SEQ_LENGTH + 2;
			maxPairEndedDistance = maxPairEndedDistance - SEQ_LENGTH + 2;
			if (pairedEndDiscordantMode) {
				maxPairEndedDiscordantDistance = maxPairEndedDiscordantDistance
						- SEQ_LENGTH + 2;
				minPairEndedDiscordantDistance = minPairEndedDiscordantDistance
						- SEQ_LENGTH + 2;
			}

			sprintf(fname1, "__%s__1", mappingOutput);
			sprintf(fname2, "__%s__2", mappingOutput);
			sprintf(fname3, "__%s__disc", mappingOutput);
			sprintf(fname4, "__%s__oea1", mappingOutput);
			sprintf(fname5, "__%s__oea2", mappingOutput);
			unlink(fname1);
			unlink(fname2);
			unlink(fname3);
			unlink(fname4);
			unlink(fname5);
		}

		sprintf(outputFileName, "%s%s", mappingOutputPath, mappingOutput);
		// Preparing output
		initOutput(outputFileName, outCompressed);

		fprintf(stdout,
				"-----------------------------------------------------------------------------------------------------------\n");
		fprintf(stdout, "| %15s | %15s | %15s | %15s | %15s %15s |\n",
				"Seq. Name", "Loading Time", "Mapping Time", "Memory Usage(M)",
				"Total Mappings", "Mapped reads");
		fprintf(stdout,
				"-----------------------------------------------------------------------------------------------------------\n");

		/********************************
		 * Regular Mode
		 ********************************/
		if (!pairedEndMode) {
			initLookUpTable();
			if (bestMode)
				initBestMapping(seqListSize);

			if (!initLoadingHashTable(fileName[1])) {
				return 1;
			}

			mappingTime = 0;
			loadingTime = 0;
			prevGen[0] = '\0';
			flag = 1;

			do {

				flag = loadHashTable(&tmpTime, errThreshold); // Reading a fragment
				curGen = getRefGenomeName();

				// First Time
				if (flag && prevGen[0] == '\0') {
					sprintf(prevGen, "%s", curGen);
				}

				if (!flag || strcmp(prevGen, curGen) != 0) {

					fprintf(stdout,
							"| %15s | %15.2f | %15.2f | %15.2f | %15lld %15lld |\n",
							prevGen, loadingTime, mappingTime, maxMem,
							mappingCnt, mappedSeqCnt);
					fflush(stdout);

					totalMappingTime += mappingTime;
					totalLoadingTime += loadingTime;

					loadingTime = 0;
					mappingTime = 0;
					maxMem = 0;

					if (!flag) {
						break;
					}
				} else if (progressRep && mappingTime != 0) {
					fprintf(stdout,
							"| %15s | %15.2f | %15.2f | %15.2f | %15lld %15lld |\n",
							prevGen, loadingTime, mappingTime, maxMem,
							mappingCnt, mappedSeqCnt);
					fflush(stdout);
				}

				sprintf(prevGen, "%s", curGen);

				loadingTime += tmpTime;
				//						lstartTime = getTime();

				initFAST(seqList, seqListSize, samplingLocs, samplingLocsSize,
						fileName[0]);

				lstartTime = getTime();

				mapAllSingleEndSeq();

				mappingTime += getTime() - lstartTime;
				if (maxMem < getMemUsage()) {
					maxMem = getMemUsage();
				}
			} while (flag);

			if (bestMode)
				finalizeBestSingleMapping();
			printf("loc_num: %lld, VHT_num: %lld \n", loc_Num, VHT_Num);
			finalizeFAST();
			finalizeLoadingHashTable();

		}
		// Pairedend Mapping Mode
		else {
			initLookUpTable();
			if (pairedEndMode)
				initBestMapping(seqListSize);

			if (!initLoadingHashTable(fileName[1])) {
				return 1;
			}
			mappingTime = 0;
			loadingTime = 0;
			prevGen[0] = '\0';
			flag = 1;

			do {

				flag = loadHashTable(&tmpTime, errThreshold); // Reading a fragment
				curGen = getRefGenomeName();

				// First Time
				if (flag && prevGen[0] == '\0') {
					sprintf(prevGen, "%s", curGen);
				}

				if (!flag || strcmp(prevGen, curGen) != 0) {

					// DISCORDANT
					lstartTime = getTime();
					outputPairedEnd();
					mappingTime += getTime() - lstartTime;
					//DISCORDANT

					fprintf(stdout,
							"| %15s | %15.2f | %15.2f | %15.2f | %15lld %15lld |\n",
							prevGen, loadingTime, mappingTime, maxMem,
							mappingCnt, mappedSeqCnt);
					fflush(stdout);

					totalMappingTime += mappingTime;
					totalLoadingTime += loadingTime;

					loadingTime = 0;
					mappingTime = 0;
					maxMem = 0;

					if (!flag) {
						break;
					}
				} else if (progressRep && mappingTime != 0) {
					fprintf(stdout,
							"| %15s | %15.2f | %15.2f | %15.2f | %15lld %15lld |\n",
							prevGen, loadingTime, mappingTime, maxMem,
							mappingCnt, mappedSeqCnt);
					fflush(stdout);
				}

				sprintf(prevGen, "%s", curGen);

				loadingTime += tmpTime;
				lstartTime = getTime();
				initFAST(seqList, seqListSize, samplingLocs, samplingLocsSize,
						fileName[0]);
				mapPairedEndSeq();

				mappingTime += getTime() - lstartTime;
				if (maxMem < getMemUsage()) {
					maxMem = getMemUsage();
				}

			} while (flag);

			if (pairedEndMode) {
				sprintf(outputFileName, "%s%s", mappingOutputPath,
						mappingOutput);
				finalizeOEAReads(outputFileName);
				outputAllTransChromosomal(transChromosomal);
				finalizeBestConcordantDiscordant();
			}

			finalizeLoadingHashTable();

			if (pairedEndDiscordantMode) {
				lstartTime = getTime();
				outputPairedEndDiscPP();
				ppTime = getTime() - lstartTime;
			}
			finalizeFAST();
		} //else

		finalizeOutput();

		fprintf(stdout,
				"-----------------------------------------------------------------------------------------------------------\n");
		fprintf(stdout, "%19s%16.2f%18.2f\n\n", "Total:", totalLoadingTime,
				totalMappingTime);
		if (pairedEndDiscordantMode)
			fprintf(stdout, "Post Processing Time: %18.2f \n", ppTime);
		fprintf(stdout, "%-30s%10.2f\n", "Total Time:",
				totalMappingTime + totalLoadingTime);
		fprintf(stdout, "%-30s%10d\n", "Total No. of Reads:", seqListSize);
		fprintf(stdout, "%-30s%10lld\n", "Total No. of Mappings:", mappingCnt);
		fprintf(stdout, "%-30s%10.0f\n\n", "Avg No. of locations verified:",
				ceil((float) verificationCnt / seqListSize));

		int cof = (pairedEndMode) ? 2 : 1;

		if (progressRep && maxHits != 0) {
			int frequency[maxHits + 1];
			int i;
			for (i = 0; i <= maxHits; i++) {
				frequency[i] = 0;
			}

			for (fc = 0; fc < seqListSize; fc++) {
				frequency[(int) (*(seqList[fc * cof].hits))]++;
			}
			frequency[maxHits] = completedSeqCnt;
			for (i = 0; i <= maxHits; i++) {
				fprintf(stdout, "%-30s%10d%10d%10.2f%%\n", "Reads Mapped to ",
						i, frequency[i],
						100 * (float) frequency[i] / (float) seqListSize);
			}
		}

		finalizeReads(unmappedOutput);
		freeMem(prevGen, CONTIG_NAME_SIZE);
	}

	return 0;
}
